Ground state structures of superparamagnetic 2D dusty plasma crystals 
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Ground state structures of finite, cylindrically confined two-dimensional Yukawa systems com- 
posed of charged superparamagnetic dust grains in an external magnetic field are investigated nu- 
merically, using molecular dynamic simulations and lattice summation methods. The ground state 
configuration of the system is identified using, as an approximation, the experimentally obtained 
shape of the horizontal confinement potential in a classical single layer dusty plasma experiment 
with non-magnetic grains. Results are presented for the dependence of the number density and lat- 
tice parameters of the dust layer on (1) the ratio of the magnetic dipole-dipole force to electrostatic 
force between the grains and (2) the orientation of the grain magnetic moment with respect to the 
layer. 
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I. INTRODUCTION 

Plasma crystals composed of dust grains that are su- 
perparamagnetic, where each grain can acquire a strong 
magnetic dipole moment in a magnetic field, are expected 
to lead to new possibilities in dusty plasma research 
[l[ . While the electrostatic interaction between the nega- 
tively charged grains is repulsive and isotropic, the mag- 
netic dipole-dipole interaction is in general anisotropic 
and can be attractive or repulsive as a function of orien- 
tation of the magnetic dipole moments 0. 

The use of superparamagnetic grains could enable 
the magnetic tuning of plasma crystal structures, sim- 
ilar to what has been considered for superparamagnetic 
colloidal crystals (e.g. Very recently, colloidal 

suspensions of sub-micron sized (« 100 nm) polyacry- 
late capped superparamagnetic magnetite (Fe304) par- 
ticles were successfully used to produce colloidal pho- 
tonic crystals with magnetically tunable stop bands cov- 
ering the visible spectrum 0t3- The superparamag- 
netic colloids form chain-like structures along an exter- 
nal magnetic field with regular inter-particle spacing, en- 
abling the diffraction of visible light. The tuning of the 
diffraction wavelength was accomplished by varying the 
inter-particle spacing. In turn this was done by varying 
the magnetic field that alters the strength of the mag- 
netic dipole-dipole interaction, which balances the repul- 
sive electrostatic interaction between the charged colloids 
0, II] • Single layer experiments with superparamagnetic 
particles on the water-air interface have demonstrated 
the advantages of the tunable inter-particle interaction 
in the studies of fundamental collective phenomena, like 
the solid-liquid phase transition [lO, Coagulation 
of charged, charged-magnetic, and magnetic dust aggre- 
gates formed from a ferrous material in various environ- 
ments was studied in (T2| , showing that the dipole-dipole 
interaction can affect the orientation and structural for- 
mation of aggregates as they collide and stick. 

While colloidal crystals typically have inter-particle 
spacings in the sub-micron regime, the spacing between 



dust grains in plasma crystals is typically larger, on the 
order of 100 fim, which is in the range of terahertz (THz) 
wavelengths. We investigate the possibility of using su- 
perparamagnetic particles in the larger micrometer size 
range in a dusty plasma monolayer in a magnetic field, 
with the aim of producing a tunable two-dimensional 
(2D) lattice structure with spacings that correspond to 
the THz regime . The tuning is accomplished by vary- 
ing the angle the magnetic field subtends with the plane 
of grains. If such structures can be produced with the 
grains occupying a large volume fraction of dust grains 
(see [l3|). they may have photonic applications in the 
THz frequency range which is currently of great interest 
owing to potential applications in spectroscopy, imaging, 
etc. [il. 

The paper is organized as follows. Section II presents 
the model system, which is a confined 2D layer of charged 
superparamagnetic grains in a plasma, placed in an ex- 
ternal magnetic field. Section HI presents the results of 
MD simulations of the ground state structures of a finite 
2D system in the crystalline solid phase as the relative 
strength of the electrostatic to magnetic dipole-dipole in- 
teraction and the direction of the grains' magnetic mo- 
ment with respect to the layer plane are varied. Section 
IV presents lattice summation calculations of the corre- 
sponding infinite 2D lattice limit of this system. A dis- 
cussion of possible experimental parameters is given in 
Section V, and a brief summary is given in section VI. 



II. MODEL SYSTEM 

The model system comprises a 2D lattice of superpara- 
magnetic dust grains immersed in a plasma in a constant, 
homogeneous external magnetic field B. Each grain ac- 
quires an electric charge q due to plasma collection, and a 
magnetic dipole moment M, which is induced by the ex- 
ternal magnetic field and therefore lies in the direction of 
B. The lattice lies in the x — y plane with an unspecified 
orientation of its principal axes. The lattice structure is 
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characterized by the lattice spacing b, the rhombic angle 
(j) and the aspect ratio v — c/b > 1, and the direction 
of its principal axes with respect to the projection of the 
magnetic field onto the plane, as shown in Fig. [T] 
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FIG. 1. Geometry of the model system, (a) The lattice lies 
in the x — y plane and its structure is characterized by the 
lattice spacing b, the rhombic angle <j> and the aspect ratio 
V = c/b. (b) The magnetic moment M of each grain lies in 
the X — z plane at an angle a to the a;-axis. The principal 
lattice axes subtends an angle /3 with the projection of the 
magnetic moment. 

The grains interact via an electrostatic screened 
Coulomb (Debye-Hiickel or Yukawa) force and by the 
induced magnetic dipole-dipole force. The electrostatic 
interaction energy between two grains with charge q sep- 
arated by a distance r is 



Ue = 9V(47reo)exp(-r/Az5)/r, 



(1) 



where A^i is the plasma Debye screening length, yielding 
a repulsive force 
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where f is a unit vector in the direction of r, which is 
the vector connecting the two particles. The magnetic 
dipole-dipole force between two grains, Fm can be re- 
pulsive or attractive depending on the relative positions 
and orientations of the grains. Since it is assumed that 
the magnetic dipole moments of all the grains are parallel 
and have the same magnitude, the interaction energy of 
the two magnetic dipoles is given by 
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The magnetic dipole-dipole force between the two grains 
is 



'M = 



47r 



[-f(5cos^6i-l) + 2mcos6'], (4) 



where f and rii are unit vectors in the direction of r and 
M, respectively, and 9 is the angle between f and m. 
In the following we choose without loss of generality our 
coordinate system such that M is in the x — z plane and 
is oriented at an angle a with respect to the a;-axis (see 
Fig.[T]). 

In a typical 2D dusty plasma laboratory experiment, 
the dust grains are confined by an electrostatic poten- 
tial. In order to approximate experimental conditions in 
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FIG. 2. (color online) Experimental (a) and MD simulation 
(b) results for the radial dust density distribution in 2D layer 
of non-magnetic dust. Distances are normalized to the ob- 
served average lattice spacing (6). The total particle numbers 
are ~ 3000 in the experiment and 5000 in the simulation. 



our simulations, we considered a more accurate repre- 
sentation of the radial dependence of the horizontal con- 
finement potential beyond the usual quadratic approx- 
imation. This was done by performing an experiment 
using non-magnetic melamine-formaldehyde (MF) parti- 
cles with the aim of measuring the radial density profile of 
the single layer dust cloud. Without going in the details, 
the experiments used 4.36 /zm diameter MF spheres, in a 
1 Pa argon gas discharge driven by 5 W of RF power at 
13.56 MHz. A layer of ^ 3000 MF spheres was created 
over the 18 cm diameter lower powered electrode. Par- 
ticle detection was performed using 650 nm wavelength 
laser illumination from the side and a 1.4 MPixel CCD 
camera from the top. Sub-pixel resolution was achieved 
using the center-of-mass method discussed in detail in 
. Calculating the density of the dust layer by averag- 
ing over the nearest-neighbor distances, the experimental 
density profile was approximated by the functional form: 



n(r) « n4r'' -I- n2r^ + no. 



(5) 



where f — r/{b) and (b) is the average lattice spacing, as 
shown in Fig.[5Ka). The corresponding molecular dynam- 
ics (MD) simulation results using N = 5000 particles and 
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an average Coulomb coupling parameter T 1000 is dis- 
played in Fig. [^b) . In contrast to infinite, homogeneous 
systems, the Coulomb coupling parameter does not fully 
characterizes the entire particle ensemble as it depends 
on the particle density, which has a strong radial profile 
in our confined system. Details of the simulation model 
can be found in the next section. 

The experimental dust density profile could be best 
reproduced assuming a horizontal confinement potential 
of the form 

V{r)^Vy + V2r^, (6) 

with V4 ^ 5-10-^ and V2 = 0.004. Here and in the follow- 
ing, we use distances normalized to the Debye screening 
length Ad, kept constant for all simulations. 

Note that compared to a confinement potential with 
a simple quadratic dependence on r, ^ this simulation 
leads to a more homogeneous distribution in the center 
of the cloud, with about 10 to 20 % lower density and 
has resulted in (b) « Ad and a dimensionless screening 
parameter k ~ (Ad V""'^) ~ 0.53 in the central region. 

III. MOLECULAR DYNAMICS SIMULATIONS 

The molecular dynamics (MD) simulations are based 
on a standard method described in e.g. [ll| . We consider 
a 2D layer particle ensemble of 5000 particles. Pair in- 
teractions (forces) are evaluated in every time-step for 
each pair of particles. Time integration is performed 
using the velocity- Verlet scheme. Particles are released 
from random positions, a slow velocity back-scaling ther- 
mostat is applied until the system reached an average 
Coulomb coupling parameter of 1000. Simulations were 
run for about 1000 plasma oscillation cycles without fur- 
ther thermostation assuming that a near to ground state 
configuration could develop during this time. 

For the simulations and the presentation of our results 
we use the following reduced quantities: Ad = 1 (length 
unit), h is the lattice spacing in units of Ad, 9 = 1 (charge 
unit) is the dust grain charge, and 77 — ixqEqM / qXo is a 
measure of the relative strength of the magnetic dipole- 
dipole interaction to the electrostatic interaction. The 
layer structure is further characterized by the bulk den- 
sity n, both b and n being an average over the central 
region of the layer, along with the rhombic angle (j) and 
aspect ratio v. 

Setting FM(r,M,a) -|-FD(r) — 0, yields the equilib- 
rium distance req(M, a) where the electrostatic and mag- 
netic forces balance. The effect of the magnetic field is 
the strongest when r is purely in the x-direction, in this 
case Eq. (g]) yields 

Fm{x, M, a) = [3(cos2a) - l] . (7) 

47r X 

In this case, there is a threshold angle, ath = 
cos~^(l/\/3) ~ 54.74°, below which the attractive inter- 
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FIG. 3. Equilibrium distance req(M,a) for a < ath along 
the a;-direction versus -q. For strong magnetic interactions 
(large rj) the attraction dominates at all distances, thus no 
equilibrium distance can be found, as indicated by the dis- 
continuation of the lines for a — 0° and 30°. 



action due to the magnetic dipole-dipole force can over- 
come the repulsive electrostatic interaction for certain 
values of 77, and agglomeration can set in. The variation 
of Toq with 77 is shown in Fig.|3]for several values of a that 
are below the threshold angle. As expected, Vgq increases 
as 1] increases, with the largest increase for small a since 
the magnitude of the attractive interaction gets larger 
as a gets smaller. Thus the grains could agglomerate at 
progressively smaller values of ry as a decreases. Further- 
more at large enough 77 values the magnetic attraction 
fully dominates over the electrostatic repulsion, thus an 
equilibrium distance can not be found at all. This trend 
will also be apparent in the following discussions of the 
MD simulation results on the variation of the structure 
of the lattice under variation of ry and a. 

An illustration of the effect of the competing magnetic 
and electrostatic interactions is shown in Fig. 21 display- 
ing the total pair potential energy U{r) = t/£;(r)-t-[/M(r) 
of a single particle for selected a angles above and below 
the threshold value and rj = 0.5. The interaction is re- 
pulsive in all directions for a — 60°; for a = 50° around 
X = an attractive region develops in a narrow angle 
around the ±x direction, separated by a potential barrier 
from the outside, as it can be seen from this cross-section 
at 7/ = (front face). Particles with high enough energy 
in the tail of the thermal distribution can overcome the 
potential barrier and result in particle agglomeration af- 
ter a long enough time. 

Turning now to the lattice structure, first consider the 
case when there is no magnetic field, so that there are 
no induced magnetic moments (77 — 0). The underly- 
ing hexagonal structure is due to the isotropic repulsive 
electrostatic interaction and is characterized hy (j) — 60° 
and V — I. Due to the boundary condition imposed 
by the cylindrical symmetry of the confinement and to 
the fact that a perfect hexagonal configuration cannot 
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FIG. 4. (color online) Total pair potential energy C/(r) = 
f/E(r) + UAi{r) around a single particle situated at r = (0,0) 
at a relative strength of the magnetic interaction to the elec- 
trostatic interaction 77 = 0.5: (a) a = 60° > ath and (b) 
a = 50° < Oth- The potential energy surface is cut at ?/ = 
to display the variation of the interaction energy along the 
X-axis. The central particle is shown. 
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FIG. 5. Radial density distribution in the dust layer (in units 
of 1/A|)) as a function of r (in units of Ad), for a — 90° and 
several values of 77. 
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form in a system with density gradinet, lattice frustra- 
tions result in slight fragmentation of the ground state 
structure. Next, consider the case there is an external 
magnetic field perpendicular to the layer, a = 90°. The 
lattice structure remains hexagonal, since both the mag- 
netic dipole-dipole and electrostatic interactions are re- 
pulsive and isotropic. As expected, the density decreases 
as rj increases (as q = 1 and Ad = 1 are kept constant), 
that is, the average lattice spacing increases owing to the 
increased repulsive force, as can be seen in the density 
profile results in Fig. El Figure |6] shows a snapshot of the 
system for a — 90° and 77 = 0.1. 

Consider now the more interesting cases when the di- 
rection of the induced magnetic moments is tilted with 
respect to the dust layer (a < 90°). The shortest lattice 
distance forms along the x-axis, thus the lattice forms 
with P = {) (see Fig[l|), as might be expected since the 
magnetic repulsion weakens or eventually turns purely 
attractive in that direction (depending on the value of 
a). Thus the system appears to align along that direc- 
tion. Figure [3 shows a snapshot of the system for a — 60° 
and rj = 0.8, where a crystal structure without domain 
fragmentation is formed, showing that the ordering effect 
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FIG. 6. Snapshot of layer for a = 90° and ?7 = 0.1 

arising from the magnetic enhancement of the interaction 
overcomes the frustration induced by the boundary con- 
dition. The variation of the lattice spacing &, the bulk 
density n, the rhombic angle (j) and the aspect ratio v 
are shown in Figs. [SlITTl respectively, as a function of -q 
for various values of a. For the large angles a > 70°, the 
lattice spacing increases and the density decreases as i] in- 
creases, because the magnetic dipole-dipole interaction is 
repulsive and its anisotropy is not strong. For this range 
of angles a, both the rhombic angle and aspect ratio of 
the lattice increase somewhat as a decreases, as the mag- 
netic interaction becomes less repulsive. For the small 
angles, a — 45° and 50°, the lattice spacing decreases 
as r] increases, owing presumably to the dominance of 
the attractive magnetic interaction, which significantly 
weakens the electrostatic repulsion. At this low a angles 
the system becomes unstable against aggregation, in the 



5 



40 I 1 , H . . J ' ■ — I 

30 - . /.■.vX:-:-:':^^'^^^^ - 
-10 

-20 -^■^v^^:^^^H^^^P^^^ 

-30 - '■•V:^i\^^'^^-:-:v^^^ " 

_4Q I 1 ' I ■ ■ • • ■ r 1 

-40 -30 -20 -10 10 20 30 40 

X / Xq 

FIG. 7. Snapshot of layer for a = 60° and tj = 0.8 



sense discussed above, at intermediate rj values. In the 
true (r — 0) ground state, low a configurations are sta- 
ble as long as 6 > , however our simulations are run at 
very low, but finite temperatures, where agglomeration 
can start (causing the simulation to stop) due to thermal 
energy fluctuations at long enough times. The rj values 
at which this occurs (77 « 0.3 ) are for this particular set 
of simulation parameters, simulation time in particular. 
During this time the system reaches only kinetically sta- 
ble states, and not a thermodynamical equilibrium state. 

This is illustrated in Fig. [T^l where the total potential 
due to the lattice particles, as experienced by a particle 
in the center is shown. In the large a > ath case (a) the 
potential energy U{r) surface exhibits a deep, well con- 
fined potential minimum. This can be contrasted with 
the case of a selected low a < ath value (b) where a min- 
imum enclosed by a low potential barrier is formed along 
the ±x directions around the vacant particle position at 
r=(0,0). 

While the central density tends to increase somewhat 
with 77, both the rhombic angle and aspect ratio increase 
rapidly, tending toward a rectangular configuration. 

At the intermediate angle a = 60°, there appears to 
be non-monotonic characteristics of some of the lattice 
parameters. When rj is small, the trends follow those 
described previously for small a, with b decreasing, and 
(j) and v increasing, as 77 increases. However, at larger 
V ^ 0.4 — 0.5, the lattice spacing b begins to some- 
what increase, although still remaining below its value 
at a = 90°. Meanwhile the central density decreases 
significantly, which may indicate that the overall mag- 
netic repulsion starts to overcome the complex effect of 
the force anisotropy. This non-monotonic behavior may 
be consistent with the trends pointed out in Fig. [S] in 
for intermediate magnetic field values, where it was 
found that for small inter-grain distances, the total force 




FIG. 8. Average lattice spacing 6 (in units of Ao) versus 77, 
for various values of angle a. 
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FIG. 9. Dust particle central density (in units of l/Af,) versus 
r;, for various values of angle a. 



Fe + Fm from eqs. ([2]) and ^ was attractive, while at 
intermediate distances the total force was repulsive. 

Overall, it can be seen from Figs. [51ITT1 that for these 
parameters, it may be possible to tune the lattice spac- 
ing and structure by changing rj and a. The lattice spac- 
ing could be tuned by a factor of about 2, with a corre- 
sponding factor of '-^ 2 in the particle density and chang- 
ing together with other lattice parameters. The lattice 
structure could be tuned from triangular (hexagonal) to 
almost rectangular; depending on rj and a, the rhombic 
angle can vary between 60 and 80 degrees, and the aspect 
ratio between 1 and 2. 



IV. INFINITE LATTICE 

We have investigated the ground state energy at T = 
(where T is the thermal energy of the dust grains) for an 
infinite, isotropic lattice when a = 90°, by using the lat- 
tice summation technique. This provides a reliable basis 
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FIG. 10. Rhombic angle versus ry, for various values of angle 
a. 
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FIG. 11. Aspect ratio v versus r;, for various values of angle 
a. 




n 


= 0.3 


n 


= 50° 




= 1 


b 


= o.es 


1/ 


= 1.45 




= 69° 



FIG. 12. (color online) Total potential energy U{r) = (7b (r) + 
17a/ (r) of the lattice experienced by a test particle situated at 
r = (0,0). Lattice parameters are taken from the MD simu- 
lations for (a) an a > ath and (b) an a < ath configuration. 
The area plotted is restricted to j/ > 0. Two particles are sit- 
uated at the front left and right corners, two other particles 
are shown at their lattice positions with their full symmetric 
(magnetic) dipole -I- (electrostatic) monopole fields. 



with high accuracy for validation of our MD simulations. 
The lattice summation was performed by summing the 
contribution of about 10® neighboring grains on a perfect 
lattice characterized by the lattice spacing 6, the aspect 
ratio v and the rhombic angle 0. In addition, in contrast 
to our previous studies of structural phase transitions 
in 2D complex plasma composed of ferromagnetic grains 
with intrinsic magnetic dipole moments .17] where the 
density was kept constant, in this case the pressure is 
kept constant as rj is varied. The pressure was computed 
from the diagonal elements of the pressure tensor, which 
in this case has the form: 

P^-lH S VrC/(r), (8) 

where 7 denotes the Cartesian coordinates {x or y), r 
is the distance between the particle at the origin (0,0) 
and another lattice particle. Summation is performed 
for particles located on a half-plane. U{r) is the inter- 
particle pair potential energy, including electrostatic and 



magnetic contributions. In the calculations the lattice is 
oriented along the x-axes, but due to the force isotropy of 
a perfect hexagonal lattice, the calculated pressure value 
does not depend on the lattice orientation. More graphi- 
cally this is the force per unit length acting on a fictitious 
line of particles inserted at x = (or y — 0) interacting 
with all particles on one side only. 

The configuration with the minimal total energy was 
sought. The initial lattice parameters were adopted from 
the = MD simulation {\d = 1, v = I, (j) = 60°, 
b = 0.81) and the initial pressure value, which was kept 
constant during the subsequent 77 > calculation, was 
evaluated for this initial lattice. The results are shown 
in Fig. 1131 which compares the lattice spacing as a func- 
tion of 77 for the infinite lattice with the MD simulation 
results for the finite system. Note that in the finite case, 
the lattice spacing is an average over the central part of 
the particle cloud. The comparison shows good agree- 
ment for the lower magnetization (jj < 0.4) cases, where 
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the deformation of the finite dust cloud is not too large 
{Ab/b ==10%) and confinement can be approximated by 
the constant pressure condition. 




FIG. 13. Lattice spacing b (in units of Xn) from "infinite" 
lattice summation and from finite MD simulations versus ry, 
for Q = 90° {f = 1, (j} = QO°). 



V. POSSIBLE EXPERIMENTAL PARAMETERS 

To aid in the design of a possible experimental realiza- 
tion of ground state structures studied in this paper, we 
estimate a range of possible plasma and dust parameters 
necessary to observe such effects. The quantity 77 is a 
figure-of-merit, which characterizes the relative strength 
of the magnetic dipole-dipole to electrostatic interaction 
between neighboring grains. Assuming that the grain is 
spherical, with radius R, and that it can be character- 
ized by a magnetic permeability /i, its induced magnetic 
dipole moment can be expressed as 



M 



B. 



(9) 



Expressing the magnitude of the grain charge as q — 
R\(l)s\ where (ps is the grain surface potential, we have 
that 



'7 = VMo^o 



M 
q>^D 



0.03 



(10) 



For example, consider a plasma with Te ~ 2 eV and 
Ui ^ 10^ cm~'^ so that the effective Debye length in the 
sheath, given approximately by the ion Debye length with 
Ti ~ Te, is about A^i ~ 1 mm. Assuming that /i = 4, 
R^5 urn, B = 5000 G, (f>s = 2 V, we obtain 77 - 0.94. 
Thus for these dust and plasma parameters, varying the 



external magnetic field from to 5000 G can vary 77 from 
to about 1. Another possibility is a denser plasma, with 
Tg ^ 2 eV and rii ^ 10^° cm"'^, and the other parameters 
the same as in the last example. In this case, varying 
the magnetic field from to 500 G can vary rj from 
to about 1. Thus it seems that there could be a range 
of reasonable experimental parameters for observing the 
variation of lattice parameters and structures predicted 
in this paper. 



VI. SUMMARY 

The ground state configuration of a 2D dusty plasma 
crystal composed of super-paramagnetic grains immersed 
in an external magnetic field has been investigated using 
MD simulations with parameters that may be close to 
realizable experimental conditions. Since the magnetic 
dipole moments of the grains are induced by the exter- 
nal magnetic field, the dipole moments of the grains all 
lie in the same direction. This study determined the de- 
pendence of the lattice parameters and structure on the 
parameter rj (which characterizes the relative strength of 
the magnetic dipole-dipole to electrostatic interactions) 
and a (the angle between the direction of the magnetic 
dipole moment and the lattice plane). It was found that, 
for a given set of dust and plasma parameters, it may 
be possible to vary the lattice spacing within a factor 
of about 2 by changing the magnitude of the external 
magnetic field or the direction of the field with respect 
to the dust layer. Gorrespondingly, the particle density 
can be varied by about a factor of 2. Moreover, the lat- 
tice structure can be tuned from triangular (hexagonal) 
to almost rectangular; depending on 77 and a, the rhom- 
bic angle can vary between 60 and 80 degrees, and the 
aspect ratio between 1 and 2. It was shown that there 
could be sets of reasonable experimental parameters for 
observing the effects discussed in this paper. There are 
other interesting dusty plasma physics issues that could 
be studied. For example, it would be very interesting to 
sec under what conditions the attractive alignment force 
may be large enough to overcome the possible rotation of 
the dust cloud due to an ion drag force in the case where 
the ions are magnetized (see e.g. [IS]). 
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